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1. Introduction 

QCD at finite density is one of the most fascinating topics of research currently because of 
the tantalizing experimental possibility of exploring the the possible region of QCD critical point 
in the T-ft phase diagram. The challenges of the conceptual difficulties, such as the complex 
fermion determinant and getting the physics of axial anomaly correct on the lattice make it further 
attractive. A lot of progress has been made in recent years [|I|] making it thus a very active field of 
research. In fhe early years of laffice QCD, one of fhe imporfanf problem was how fo incorporafe 
a finite chemical potential, ft, in the local Dirac operator on the lattice. Naively introducing jx 
as a Lagrange multiplier of the local point-split conserved charge ^fx74Yx+6 + Wx+67^¥x led to /r- 
dependent divergences in the free energy of quarks in the continuum limit. Arguing that since 
/I enters in the continuum like a fourth component of a static gauge field, fhe chemical pofenfial 
should enter on the lattice as e^^® factors, multiplying the forward and backward temporal gauge 
links respectively 0. It did eliminate the undesired divergences mentioned above. However, this 
proposal is not unique though: weights (1 ±a/r)/y^l —a^pL^) to the temporal links also lead to 
finite results for the free energy [^]. In fact the most general weights /(/ra),g(/ia) should follow 
/{aji.) ■ g{aiJL) = 1 with /(O) = f'{0) = 1 in order to cure the divergences on the lattice [Q]. The 
analytical proof of absence of divergences in all these cases above was for free quarks. Indeed, 
further numerical computations in quenched QCD showed it to work for the interacting case as 
well lH], while a similar check in the full theory is still lacking. 

We revisit the problem of introduction of chemical potential on the lattice again with a view 
to develop a simpler alternative. We foresee two important applications in which the alternative we 
pursue here may be of practical relevance in the coming years. Heavy ion collision experiments at 
RHIC at Brookhaven National Laboratory and at FAIR, GSI may likely produce strongly interact¬ 
ing systems at finite density. An important theoretical input for modeling of the evolution of these 
systems is the QCD equation of state(EoS). In order to calculate the EoS as a function of one 
needs higher order quark number susceptibilities. The latter are also important measurables on the 
lattice which relate directly to the fluctuations of net baryon number that are being measured 
in the heavy ion collision experiments. Eurthermore, these are also important in estimating the 
location of the critical point in the QCD phase diagram. The critical point is governed by the sin¬ 
gularity of the Taylor expansion of the baryon number susceptibility in Pb/T. Its location can be 
determined by the radius of convergence of this Taylor expansion. It, in turn, depends on the ratios 
of the higher order quark number susceptibilities. Therefore, their accurate and efficient measure¬ 
ment on the lattice is necessary although a rather challenging problem. Having a Dirac operator 
with a linear pN term, simplifies fhe problem quife a lof. Eor example, fhe number of fhe inver¬ 
sions of fhe Dirac operafor needed fo calculafe fhe eighfh order quark number suscepfibilify wifh 
fhe naive sfaggered quarks reduces from 20 fo 8 0. Eor ofher improved Dirac operafors and even 
higher order QNS, fhe gain is even more which leads fo a significanl reducfion of fhe compufafional 
time. 

Since fhe sfudy of QCD crifical poinf would be a major problem of inferesf in fhe coming 
years bofh from experimenfal and fheory perspecfive, if is also imporfanf fo consider how crucial 
is fhe role of exacf chiral fermions af finite densify for fhe laffice sfudies. The presence of fhe 
crifical end-poinf depends crucially on fhe number of lighf quark flavours. Model-based consid- 
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erations favour its presence for QCD with only two light quarks [||]. In addition, it needs the 
anomalous f/A(l) symmetry-restoration to take place at sufficiently high temperatures [|^. A mod¬ 
erately heavy strange quark may affect the location of the QCD critical end-point quantitatively but 
not its existence. Majority of the calculations on QCD at finite density on the lattice employ stag¬ 
gered fermions and its improved versions. However, the continuum flavour and spin symmefries 
are infermingled for the staggered quarks and the flavour single! Ca(1) anomaly is recovered only 
in the continuum limit. As is well-known, the overlap and/or domain wall fermions [ [I^ are 
much more preferable from the chiral symmetry perspective. They have both the correct chiral and 
flavour symmefry on lattice as well as an index theorem on the lattice |j^ [^. These are likely 
to be crucial for investigations of the QCD critical point. At present the computations with chiral 
fermions are prohibitively expensive. 

However, recently domain wall fermions have been used to measure the chiral crossover tran¬ 
sition temperature on the lattice [|I^. This gives us a hope that with the increase in computational 
resources and smarter algorithms, use of chiral fermions would be more realistic in the coming 
years. Non-locality of the overlap fermions makes the introduction of the chemical potential non¬ 
trivial. Bloch and Wettig [M] proposed to use the same exponential prescription as above for the 
timelike links of the Wilson-Dirac kernel Dwiaft) to define fhe corresponding overlap Dirac ma- 
frix af nonzero density. The free energy from fhis overlap fermion acfion has no divergences 
lElllli in fhe free case. Unfortunalely, however, if has no chiral invariance for nonzero /r eifher 
[jl^. Using fhe definition of fhe chiral projecfors for overlap fermions, we |[T^] proposed a chirally 
invaiianf Overlap acfion for nonzero jJ. : 


n 

— [^^<9V ~l~ 1 • 

n 

It was shown that the fermion action is invariant under the Luescher transformation 0 
and fhe corresponding order paramefer, fhe chiral condensafe in unique for all values of afx [JT^]. 
Moreover fhe anomaly is /r-independenf on fhe lattice as is expecfed from fhe continuum [|T8|]. If, 
however, has fhe same /i-dependenf divergences in fhe number density and fhe energy density 
as fhe linear /r-case for naive/sfaggered fermions JT8[]. Furthermore, unlike fhaf case, fhese cannof 
be removed by exponenfiafion of fhe /i-ferm [jl^. So eifher one has exacf chiral invariance on 
fhe lattice or have fo deal wifh fhe divergences in fhe confinuum limif of a —> 0. One needs fo 
undersfand fhe origin of fhese /r-dependenf divergences beffer. In quanfum field fheories, removal 
of divergences from physical quanfum theories has lead to extensive studies of the different reg¬ 
ulators and their suitability for different physical problems. We aim to have an understanding of 
the nature of the divergences at finite density, and look for methods to remove them in alternate, 
perhaps simplifying, ways. 


2. The divergences at finite density : Why and how to remove them 

2.1 Non-interacting fermions at finite chemical potential 

As a first step we examined carefully the free dense quark gas in continuum. We found that 
contrary to the common belief, these divergences are not due to lattice artifacts. Indeed, the jx- 
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dependent divergences exist in the continuum theory as well when appropriate care is taken in 
manipulating divergent integrals. The lattice regulator simply makes it easy to spot them. Using 
a Pauli-Villars cut-off A in the continuum theory, one can also show the presence of terms in 
the number density easily [|20[]. We will sketch the argument here briefly, refening the reader to 


[ ]20| ] for more details. Let us recall that QCD partition function can be written in the path integral 
formalism as, 

= J (2.1) 


where (y) V and represents the (anti-)fermion and the gluon fields respectively. One introduces 
a chemical potential /r, for each conserved charge Alf = f d^xYiY4Yi’ leading to the additional term 
in the action above. Canonical definitions consisting of first and second derivative of ^ with 
respect io iiq = ft yield quark number density and susceptibility. 


T d\n^{ft). 


V 


T—fixed 


( 2 . 2 ) 


These can be calculated analytically when the gauge interactions are switched off. For simplicity 
we consider only massless fermions though this derivation can be easily extended to finite quark 
mass. The expression for the free quark number density is then 



d^p {(On + ilJ.) 
ilTlf p'^ + {(On + i}xY 



d^p 

(271)3 


Y,F{(On,P,p), 


(2.3) 


where p^ = p^+P 2 +pI- All the gamma matrices are Hermitian in our convention. T = 0 and /i = 0 
corresponds to the vacuum contribution which can be removed by subtracting n{T = 0,jU = 0). 
Although this is identically zero for number density, the corresponding subtraction for the energy 
density is actually oc A^. Due to the Fermi-Dirac distribution functions, one does not expect any 
ultraviolet divergences at finite T. We therefore consider the T = 0 contribution of the above 
expression to examine the presence of divergences, if any. The sum over &)„ turns into a continuous 
ft) integral and the explicit expression for the zero-T contribution to the number density is. 


/“^ d(0 d^p {(O + ip.) 

” J -00 2n {2nY p^ + {(0 + ipY ' 


(2.4) 


Under a variable transformation co + ip = ft)', it can be recast as 


n = 4i 


dco' d^p 


(O' 


(2.5) 


2 n (271)3 _|_ ^/2 • 

In calculating this expression one considers a suitable contour as shown in Figure with a cut-off 
A for the (O integral in order to regulate the divergent terms, and identify them by powers of A. The 
expression in Eq. ^ can be seen to be the line integral 3 in Figure ||. By adding and subtracting 
the line integrals labeled as 1,2,4, one obtains the result from the contour integral in terms of the 
residue of the quark propagator pole with the corresponding 0-function defining the Fermi surface 
and leading to the usual /i3-term. Out of the subtracted terms, the line integral 1 is exactly zero 
due to CP symmetry. However, the contribution of the arms 2 and 4 in the Figure clearly do not 
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Figure 1: The contour integral used for calculating the zero temperature part of the number density. 


cancel each other. The terms can be shown to arise from the algebraic sum of these terms [p^], 
being of the form log . In the limit jj., and with the same cut-off A for momentum 

p, an expansion in p/A shows that while the leading A^ terms do cancel from the numerator and 
denominator, the pA^ terms add, and survive. It is easy to see that a similar p^A^ divergent term 
persists in the energy density as well. It is amusing to note that the divergence is absent in the 
number density for nonzero isospin chemical potential with jUg = 0, since this amounts to putting 
two flavours of quarks in the Lagrangian with chemical potentials jU/ and —/f/. The divergent term 
then cancels explicitly between the two quark flavours since it comes with an equal and opposite 
sign. The corresponding energy density, however, does have the p^A^ divergent term. On the other 
hand, the divergence does exist for imaginary chemical potential for the number density as well, as 
we show below. 


2.2 Non-interacting fermions in presence of imaginary p 

The expression for the number density for free fermions in presence of imaginary chemical 
potential at zero temperature is 


d(0 d^p (ft) — p) 

-A 2n (271)3 p^ + {(o- pY ' 


( 2 . 6 ) 


Performing the ft) integral, which in this case is on real ft)-axis only, one gets. 


n = —4/ 


I 


d^p 


In 


/ p2 + (A-^)2 X 

\p^ + {a+pyJ 


(2.7) 


In the limit A S> ft, one can again expand the numerator and the denominator of the integrand. As 
above, we further assume the same cut-off A for the momentum integral. 

Keeping terms up to l/A^, one obtains the number density as. 


p^ pA^ 


( 2 . 8 ) 


The imaginary nature of number density comes out, as expected. Indeed, one recovers the same quadratic 
divergence back if one perform analytic continuation to real p. This is reassuring since it shows 
that the presence of divergence for nonzero p is not due to the particular contour method we chose 
in Figure [I], as the integrals in this case are performed on the real ftj-axis for imaginary chemical 
potential. 
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2.3 Free quark gas on the lattice 

We have seen that for non-interacting fermions within a cut-off regulator scheme one obtains 
a jUA^ divergence. This explains the presence of a similar divergence in the lattice expression of 
number density for non-interacting fermions, since the lattice regulator is simply a cut-off regulator 
and the /r/a^ term in the number density on a lattice is just a manifestation of the continuum 
divergence. This can be explicitly seen following the same computation as above. The expression 
of the number density on a x Nt lattice with linear /r and a conserved point-split charge is, 

3 / y-, (sin (On + /ft a cos (O,,) cos to,, 

flQ- = — - - / - - -. (2.9) 

^ /+(sinoj,,+ ZjUacosoj„)^ 

where/= {maf-+ fiv^{ap\) + ?,m^{ap 2 ) + ?>m^{apf),pj = lnj%/N^ nj gZ. AtT = 0 one converts 
the sum over the Matsubara frequencies into an integral, and obtains a line integral in the complex 
ft) plane analogous to the line integral marked 3 in Figure ||. It can be calculated using a contour 
as in Figure but with tanh^^(a/t) in place of ft. Again the residue of the contour integral yields 
the Fermi surface on a finite lattice even in this method. The divergence now comes from the line 
integral equivalent to the one marked 1 in Figure On the lattice one can perform the integral 
over ft) to obtain the zero temperature part of the number density. Taking successive derivatives 
of this term with respect to ft gives the zero temperature artifacts in the linear ft method for the 
second and fourth order susceptibilities [ |I^ : 

( 2 . 10 ) 


3. Testing the idea on the lattice 


Having noted that the divergence on the lattice is no different that a A-cut-off in the continuum 
theory, we proposed, and demonstrated, earlier that a subtraction of the T = 0 ideal gas term for 
the quark number susceptibility on the same size lattice does eliminate the divergence in the 
continuum limit obtained by sending Nj ^ for the free case [0- Here we concern ourselves 
with the test of that idea for interacting theory. We calculated the susceptibilities(QNS) for two de¬ 
generate quark flavours on quenched lattices with the linear ft Dirac operator p0|]. The pure SU (3) 
gauge configurations were generated using the Cabibbo-Maiinaii pseudo-heatbath algorithm with 
three SU (2) subgroup update per sweep using the Kennedy-Pendleton updating method. We gen¬ 
erated Nt =4,6,8,10 and 12 lattices and 25-100 independent configurations each at two different 
temperatures given hy T/Tc = 1.25,2. To keep the finite volume effects under control we checked 
and used Ns = 4 Nt - The second order QNS was calculated by inverting the Dirac operator on 400 
and 500 Gaussian random vectors at 2,1.25rc respectively. The results obtained using this method 
contain unphysical ft^/a^ and artifacts for second and fourth order QNS. For removing them 
we calculate the corresponding free theory QNS at zero temperature in Eq. ( p.lO ) and perform the 
continuum extrapolation. The expressions for these T = 0 values and the corresponding numeri¬ 
cal values are displayed in Table ||. Performing the subtraction, the physical values of the second 
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Nj 

r(MeV) 

miTc 

X2o{T = 0) 

;f4o(r = o) 

4 



-0.062487 

-0.142923 

8 

1.25 


-0.062499 

-0.142661 

10 


0.1 

-0.062501 

-0.142620 

12 



-0.062502 

-0.142607 

4 



-0.062541 

-0.142999 

6 

2 


-0.062500 

-0.142611 

8 


0.1 

-0.062508 

-0.142677 

10 



-0.062512 

-0.142718 


Table 1: The T = 0 ideal gas subtraction terms at different volumes and temperatures for the second and 
fourth order QNS used in this work. 




Figure 2: The continuum extrapolated results for the second order quark number susceptibilities in quenched 
QCD for a quark mass m/Tc = 0.1 at 1.257’^ (left) and 2Tc (right) compared to the corresponding results using 
the e^^ method. 


order QNS with a quark mass m/Tc = 0.1 is shown in Figure ^ We do not observe any divergent 
term as evident from the positive slope of the data. Moreover, our extrapolated continuum result 
coincides with the earlier result obtained with the exp(±ajU) action [||, ^|, marked with red open 
circles on the vertical axis of each figure. In other to check that there are no mass dependent di¬ 
vergent terms we calculated the second order QNS with a still lower quark mass m/Tc = 0.01, the 
continuum extrapolated results shown in the left panel of Figure We indeed confirm absence 
of any new divergences from fhe slope of fhe filled curve. Furfhermore, removing fhe free Iheory 
r = 0 arlifacls from fhe fourlh order QNS also gives fhe correcf continuum limil as evidenl from 
fhe righl panel of Figure ^ For fhe QNS Xn ,n> 6, fhe additional arlifacls in Ihe linear jj. melhod 
are hence do nol affecl Ihe continuum exlrapolalion of Ihe lallice resulls. Moreover Ihese 

terms reduce Ihe cul-off effecls compared lo Ihe e^a-melhod and facililale a smoolher approach lo 
Ihe continuum [0^]. 


4. Summary 


We did a comprehensive sludy of Ihe possible origin of Ihe divergences in Ihe quark 
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Figure 3: The continuum extrapolated results for the second order quark number susceptibilities for a lighter 
quark mass m/Tc = 0.01 at 1.257’c (left) and the results for the fourth order susceptibility for the quark mass 
ot/Tc = 0.1 (right). 

number susceptibility on the lattice when we introduce chemical potential as a ptN term in the 
fermion action, N being the conserved number density on the lattice. We found out that such 
a divergence explicitly exist for the QNS of free fermions with a cut-off regulator. The lattice 
regulator just faithfully reproduces this divergence for the free fermions as well. Conserved charges 
are not renormalized. Assuming therefore that the only divergence that would exist for QCD on 
the lattice is due to the ptN term of the free theory, we did an explicit check in quenched QCD by 
taking the continuum limit of the second and the fourth order ideal gas subtracted QNS. We indeed 
find that there are no additional divergences that arise due to interactions. The divergence in the 
linear jU-method can be removed easily which gives us another path towards introducing jx on the 
lattice in addition to the most popular method [Q] where one modifies fhe quark acfion af finite jx 
to remove these divergences explicitly. This method may be beneficial for calculafing fhe higher 
order QNS wifh considerably less compufafional efforf allowing us fo calculate fhe Laffice EoS af 
finife baryon chemical potential and make progress fowards measuring fhe critical end-poinf. 
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